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We perform a detailed study of the weak interactions of standard model neutrinos with the 
primordial plasma and their effect on the resonant production of sterile neutrino dark matter. 

Motivated by issues in cosmological structure formation on small scales, and reported X-ray signals 
that could be due to sterile neutrino decay, we consider 7 keV-scale sterile neutrinos. Oscillation- 
driven production of such sterile neutrinos occurs at temperatures T > 100 MeV, where we study 
two significant effects of weakly charged species in the primordial plasma: (1) the redistribution of 
an input lepton asymmetry; (2) the opacity for active neutrinos. We calculate the redistribution 
analytically above and below the quark-hadron transition, and match with lattice QCD calculations 
through the transition. We estimate opacities due to tree level processes involving leptons and 
quarks above the quark-hadron transition, and the most important mesons below the transition. 

We report hnal sterile neutrino dark matter phase space densities that are signihcantly influenced 
by these effects, and yet relatively robust to remaining uncertainties in the nature of the quark- 
hadron transition. We also provide transfer functions for cosmological density fluctuations with 
cutoffs at A: ~ 10 h Mpc“^, that are relevant to galactic structure formation. 

PACS numbers: 95.35.+d,14.60.Pq,14.60.St,98.65.-r 


I. INTRODUCTION 

Deep in the radiation dominated epoch of the Universe, 
the three neutrinos present in the standard model (SM) 
of particle physics ^ make up a significant population 
of relativistic species within the primeval cosmic plasma. 
We have strong evidence of their existence at these early 
epochs from probes of the primordial Universe such as 
the Cosmic Microwave Background (CMB) (probing tem¬ 
perature T ^ 0.25 eV) [2], and the synthesis of light 
elements during the epoch of Big Bang Nucleosynthesis 
(BBN) which depends on the neutron-to-proton ratio set 
at Tdec ~ 1.5 MeV [3], the temperature of weak neutrino 
decoupling. Above this temperature, SM neutrinos inter¬ 
act with species that carry weak charge, through which 
they remain coupled to the primordial plasma [1]. 

There is a long history of speculation about additional 
neutrino species (see Ref. [5] for a recent review). Owing 
to the precise measurement of the invisible decay width 
of the SM Z boson [T] , any extra neutrino species must be 
“sterile” (i.e. electroweak singlets) [5]. Furthermore, pre¬ 
cise measurements of the CMB HE] and of the primeval 
abundance of light elements [iin] strongly constrain the 
presence of extra relativistic species in the early Universe. 
These constraints indicate that (i) unlike SM neutrinos, 
light sterile neutrinos never fully thermalize with the rest 
of the cosmic plasma unHH], or (ii) that sterile neutri¬ 
nos are massive enough to form the inferred population 
of dark matter (DM) in the Universe (see e.g. Ref. |15)). 
Sterile neutrinos with masses in the keV range act as 
DM in the CMB era, but are relativistic in the BBN era, 


when they do not significantly impact the expansion rate 
due to their negligible energy density (compared to the 
Fermi-Dirac value). 

Early works in this direction studied right-handed ster¬ 
ile neutrinos with masses rUg « 0.1 — 100 keV, produced 
by the oscillation of left-handed SM neutrinos [TBHIO]. 
The mixing angle between the SM and sterile neutrinos 
is fixed by the present day DM abundance. In the original 
Dodelson-Widrow scenario HZl, sterile neutrinos are pro¬ 
duced with a momentum distribution reflecting that of 
the active neutrino species, and thus constitute “warm” 
DM pTM24j . However, small-scale structure formation 
and X-ray observations [55H57] appear in signifi¬ 
cant conflict with the fiducial Dodelson-Widrow scenario, 
hence prompting the search for alternative sterile neu¬ 
trino production mechanisms [T51 155H55] . 

In this paper, we examine in detail the resonant pro¬ 
duction of sterile neutrinos in the presence of a small 
primordial lepton asymmetry. Originally proposed by 
Shi and Fuller |38j . this production mechanism makes 
use of a small lepton asymmetry to modify the plasma’s 
interaction with SM neutrinos in such a manner as to res¬ 
onantly produce sterile neutrinos at particular momenta 
[Him US]- This generically results in ‘colder’ DM dis¬ 
tribution which improves consistency with models of cos¬ 
mological structure formation |57Hf)5j . while requiring a 
modest primordial lepton asymmetry, which is relatively 
poorly constrained [55H7D] . 

Sterile and active neutrino mixing, which is needed for 
the former’s production, also leads to their decay m 
Eg. For typical values of the sterile neutrino mass this 
predicts an X-ray flux from the DM distribution in the 
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low redshift Universe This has been the subject 

of much recent interest, due to hints of an excess flux 
at ^ 3.5 keV in stacked X-ray spectra of several galaxy 
clusters m and in observations of M31, the Milky Way, 
and Perseus [ZSITS]. There is currently an active debate 
on the existence, significance and interpretation of this 
excess [7S1 - I55] . In the present work, we use this tentative 
signal as a motivation to study in detail the physics of 
sterile neutrino production in the early Universe, but the 
machinery we develop is more generally applicable to the 
broader parameter space of the Shi-Fuller mechanism. 

We present here an updated calculation of resonantly- 
produced sterile neutrinos and relax several simplifica¬ 
tions that had been adopted previously in the litera¬ 
ture. Furthermore, we leverage recent advances in our 
understanding of the quark-hadron transition in order 
to include a more realistic treatment of the strongly- 
interacting sector. Our motivation is twofold: a) improve 
the treatment of lepton asymmetry, which is a crucial 
beyond-SM ingredient in the mechanism, and b) provide 
realistic sterile neutrino phase space densities (PSDs) 
and transfer functions for matter fluctuations, which are 
starting points for studying cosmological implications on 
small scales. Our improvements to the sterile neutrino 
production calculation can broadly be classified in three 
categories. 

Firstly, we study how the cosmic plasma reprocesses a 
primordial lepton asymmetry. For models that can ex¬ 
plain the above X-ray excess, the majority of sterile neu¬ 
trinos are produced at temperatures above 100 MeV [S0] . 
At these temperatures, there is a significant population 
of either quarks or mesons, depending on whether the 
temperature is above or below the quark-hadron transi¬ 
tion. Since these hadronic species are coupled to neu¬ 
trinos and charged leptons through weak processes, the 
establishment of chemical equilibrium among the differ¬ 
ent constituents of the cosmic plasma will automatically 
transfer a primordial lepton asymmetry to the hadronic 
sector. An illustrative example is the reaction 

+ — 7r+, (1) 

which can redistribute an initial neutrino asymmetry into 
charged lepton and hadronic asymmetries. At lower tem¬ 
peratures, the asymmetry is redistributed to a lesser de¬ 
gree between the leptonic flavors. As we discuss in the 
body of the paper, this redistribution modifies the dy¬ 
namics of the resonant sterile neutrino production, re¬ 
sulting in a modified final PSD. 

Secondly, we incorporate several new elements to the 
calculations of the neutrino opacity (i.e. the imaginary 
part of the self-energy) at temperatures 10 MeV < T < 
10 GeV. Accurate neutrino opacities are needed since 
they basically control the production rate of sterile neu¬ 
trinos through cosmic epochs. Early works on neutrino 
interactions in the early Universe [4l[T8l|84] assumed that 
neutrinos largely scattered off relativistic particles and 
thus scaled their cross-sections with the center-of-mass 
(CM) energy. In addition, these calculations also ne- 
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FIG. 1: Sterile neutrino DM parameter space: shaded 
regions are consistent with the X-ray signal at 1, 2 and 3 a. 
The best determined parameters are from the MOS stacked 
clnsters of Ref. m- Statistically consistent signals are found 
in their core-removed Perseus spectrum, and M31 |74| . The 
lines show constraints at the 90% level from Chandra 
observations of M31 (Hid) |87| . stacked dwarf galaxies (M14) 
[88j . and Suzaku observations of Perseus (T15) |89j . Stars 
mark the models that we study in the body of the paper. 


glected the effects of particle statistics. Under these two 
simplifying assumptions, the opacity r(Ay^) for an input 
neutrino of energy Ei,^ is of the form 

T{E,^) = \{T)GlT^E,^, (2) 

where Gp is the Fermi coupling constant, and A(T) is a 
constant that depends on the number and type of avail¬ 
able relativistic species in the cosmic plasma. References 
POI155] subsequently developed a framework to include 
particle masses, loop corrections, and particle statistics 
in the neutrino opacity calculation. In the present work, 
we add previously-neglected contributions to the opacity 
such as two- and three-body fusion reactions, and also 
use chiral perturbation theory to compute the hadronic 
contribution to the opacity below the quark-hadron tran¬ 
sition. We find both quantitative and qualitative mod¬ 
ifications to the form of Eq. (i). Wherever we present 
matrix elements, we use the ‘— I—I—1-’ metric signature. 

Thirdly, we fold the asymmetry redistribution and 
opacity calculations into the sterile neutrino production 
computation, and provide updated PSDs for the range of 
parameters relevant to the X-ray excess. As part of this 
process, we carefully review and correct the numerical 
implementation of the sterile neutrino production used in 
Ref. [53] . Our sterile neutrino production code is publicly 
available at https://github.com/ntveem/sterile-dm 
We finally use the updated sterile neutrino PSDs in a 
standard cosmological Boltzmann code [86] and provide 
new dark matter transfer functions. 

We organize the paper such that the beginning sections 
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deal with SM physics, while the later ones apply their 
results to sterile neutrino DM. We introduce the produc¬ 
tion mechanism in Se c.[IT| We then study the asymmetry 
redistribution in Sec. |III| and active neutrino opacities in 
Sec.|IV] Finally, we apply these results to sterile neutrino 
production in Sec.|Vj evaluate transfer functions for mat¬ 
ter fluctuations in Sec. lYg and finish with a discussion of 
our assumptions and uncertainties in Sec. lYTTl We collect 
technical details into the appendices. 


II. OVERVIEW OF RESONANT STERILE 
NEUTRINO PRODUCTION 

In this section, we briefly review the resonant pro¬ 
duction of sterile neutrinos in the early Universe. We 
first present the specific scenario that we consider in this 
work, and then discuss the Boltzmann formalism used 
to compute the out-of-equilibrium production of sterile 
neutrinos. We finally discuss how the presence of the 
thermal bath and lepton asymmetry change the neutrino 
self-energy and govern the sterile neutrino production. 
We refer the reader to Refs. [T8l[38l|56] for more details. 


A. Assumptions 


In our current study, we focus on the following sce¬ 
nario. 


1. We consider an extra sterile neutrino species, Vg, 
that is massive compared to the SM neutrinos 
which we take to be effectively massless. 
The propagation (light/heavy) and interaction (ac¬ 
tive/sterile) eigenstates are related by a unitary 
transformation, the most general version of which 
is a 4 X 4 matrix. We assume that the sterile neu¬ 
trino mixes with only one of the SM ones, which 
we take to be the muon neutrino, i.e., 


_ ( cos6» sin6»\ / 4'o \ 
\-sin0 cosO) J ■ 


(3) 


The fields on the left- and right-hand sides are in¬ 
teraction and mass (rus) eigenstates, respectively, 
and 9 is the active-sterile mixing angle. The choice 
of a muon neutrino is arbitrary, and reflects the 
choice of previous work [181 US] ■ 

2. We assign a non-zero lepton asymmetry to the pri¬ 
mordial plasma. In the general case, each SM flavor 
has its own asymmetry, but we assume a non-zero 
value only for the mu flavor (i.e. the one that mixes 


with the sterile neutrino): 

T — ^QL (4) 

where the dimensionless asymmetry An^ in species 
A is the temperature-scaled difference between the 
particle and anti-particle densities, ISua = (riA — 
nA)lT^, and 5afi is the Kronecker delta. In general, 
entropy-scaled asymmetries are preferable, since 
they are conserved through epochs of annihilation. 
However, the definition used in Eq. Q simpli¬ 
fies comparison with lattice QCD calculations in 
Sec. |III[ We fix by hand the mu lepton asymmetry 
at high temperatures to produce the canonical DM 
density, = 0.1188 in the current epoch [5]. 

In the rest of the paper, we use a hat to indicate tem¬ 
perature scaled quantities. We choose to study the pa¬ 
rameter space of interest for resonantly produced ster¬ 
ile neutrino DM consistent with the recent X-ray signal. 
Figurej^shows a section of the rUg and sin^ 29 plane with 
contours for the unidentified lines of Refs. [731 El], along 
with constraints from Chandra observations of M31 [57], 
stacked dwarf galaxies [55], and Suzaku observations of 
Perseus [M] . The stars show a range of mixing angles at 
a specific value of rUs , and mark models that we study in 
Sections |V| and IVII 

For all these models, the bulk of the sterile neutrinos 
are produced at temperatures well below the masses of 
the weak gauge bosons 80 GeV), but above weak de¬ 
coupling at T ~ 1.5 MeV [SS]. Active-active neutrino 
oscillations in the primordial plasma are suppressed at 
these temperatures [90], hence it is consistent to assign 
individual asymmetries in Eq. Q and neglect electron 
and tau neutrino mixing in Eq. 


B. Boltzmann Formalism 

In its full generality, out-of-equilibrium sterile neutrino 
production (via oscillations) is best described by the evo¬ 
lution of the two-state density matrix of the neutrinos in 
the active-sterile (interaction) basis pTI194| . 

For the parameter range in Fig. [^ most sterile neutri¬ 
nos are produced above temperatures T > 100 MeV. At 
these temperatures, the two state system is collisionally 
dominated, i.e. the interaction contribution dominates 
the vacuum oscillations. In this regime, the evolution 
of the density matrix separates out and yields a quasi- 
classical Boltzmann transport equation for the diagonal 
terms, which are the PSDs of the active and sterile com¬ 
ponents [HSHST]. The Boltzmann equation for the sterile 
neutrino PSD is 


m 


fvAvA) - Hp 


dp 


fuApA) = 
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E 


d^Pa 


d^Pi 


. . . {2tt) S {p+Pa-\ - P^ - ■■■) 


- 


1 

X - 

2 L 


{2TTf2Ea ■ ■ ■ (27r)32£;, 

{Pm(yfi —>■ Vs\p, t)) (1 ~ fv^) - ya+v^ + ...fi • ■ • (1 T fa) (l ~ /i^^) • ■ ■ 

— {Praia's ^^l\P^ t))fvg (l ~ fu^) |A^|i,^_|_q_|- yi+...fa • ■ • (1 T /i) ■ ■ ■ 


(5) 


We can write an analogous equation for the antineu¬ 
trinos. Here, the f{p) are PSDs for particles with three- 
momentum p and energy E, and H is the Hubble expan¬ 
sion rate. The right-hand side sums over all reactions 
that consume or produce a muon neutrino. The symbol 
^|A4p denotes the squared and spin-summed matrix 
element for the reaction, and the multiplicative factors 
of (1 =F /) implement Pauli blocking/Bose enhancement 
respectively. The factor of 1/2 accounts for the fact that 
only one (i.e. the muon neutrino) state in the two-state 
system interacts The Pm are active-sterile os¬ 

cillation probabilities in matter, which depend on the 
vacuum mixing angle 0, and are modified by interactions 
with the medium. The latter are parametrized by the 
neutrino self energy [M], and the ‘quantum damping’ 
rate for active neutrinos. In terms of these quantities, 
the oscillation probabilities are [Si [57] 

^ Vs\P,t)) = {l/2.)/\^{p)s\V?29 

X |A2(p)sin2 26l-f 

-k [A(p)cos26»-H^-H*‘’(p)]^} \ (6) 

We have introduced the symbol A(p) for the vacuum os¬ 
cillation rate, A(p) = 5m‘f,^ ^^^/2p, and split the neutrino 
self energy into the lepton asymmetry potential and 
the thermal potential (the asymmetry contribution 
enters with the opposite sign in the version of Eq. ®> 
for antineutrinos). The quantity D{p) is the quantum 
damping rate, and equals half the net interaction rate of 
active neutrinos [the factor of half enters for the same 
reason as it does in Eq. ([^]. The net interaction rate for 
a muon neutrino is 




E 






dfpa dfpi 

(27r)32P, • ■ ■ (2^)32P, ’'' 


X {2tt)^5'^{p + Pa^ - Pi - ...) 


^ E \^\l^+a+■■■^^+.Ja . . . (1 T /z) • ■ ■ 

(7) 


We simplify the phase-space integrals in Eq. ([^ by us¬ 
ing detailed balance to equate the forward and back¬ 
ward reaction rates. The resulting Boltzmann equation 
for quantum-damped and collisionally-driven sterile neu¬ 
trino production is | 18 j 


^UAP^t) - Hp -^f^^{p,t) 

~ ^ u^;p^t)) [U^{p,t) - UAP,t)] : 

( 8 ) 

with a related equation for antineutrinos. There are sub¬ 
tleties with the effects of quantum-damping in the case 
of resonance [55], but tests with the full density matrix 
formalism find that the quasi-classical treatment is ap¬ 
propriate [56]. 


C. Asymmetry and Thermal Potentials 


We now expand on the origins of the asymmetry and 
thermal potentials appearing in Eq. ([^. These potentials 
encapsulate the self energy of propagating active neutri¬ 
nos due to interactions with the plasma. Under the con¬ 
ditions we are interested in, there are three contributions 
to the neutrino self energy: a) an imaginary part pro¬ 
portional to the net neutrino opacity, b) a real part due 
to finite weak gauge boson masses (U*^), and c) a real 
part proportional to asymmetries in weakly interacting 
particles (V^). We follow the treatment in Ref. [S3], and 
recast it in terms of the quantities that we compute later. 

Figure [^ shows lowest-order contributions to active 
neutrinos’ self energy. Thick red lines are thermal prop¬ 
agators of weakly charged species in the background 
plasma. There are two corrections - bubbles and tad¬ 
poles, shown in Fig. 2a and |2b| respectively. The back¬ 


ground fermion is a lepton of the same flavor in the for¬ 
mer, and any weakly charged species in the latter. 

A massless active neutrino’s ‘dressed’ propagator is 


iPi'c ) = = -Pi'c ■ (9b) 


Here, p,^^ and u are the neutrino and plasma’s four- 
momenta, A is shorthand for and is the left 


handed neutrino’s self energy. Equation (9b) divides this 


self energy into two contributions that affect the particle 
and anti-particle poles of Eq. (9a) differently. Figure 


illustrates their association with asymmetry and thermal 
potentials: 





(10) 
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(a) bubble diagram (b) tadpole diagram 


FIG. 2: Lowest order contributions to a propagating active 
neutrino’s self energy. Red lines are thermal propagators. In 
(a), / is any species with weak charge. In (b), f = Va, a~. 



FIG. 3: Matter potentials for massless neutrinos in the 
plasma’s rest frame: filled and un-filled circles are poles at 
finite and zero temperature, respectively. See Sec. |IV| for the 
imaginary shift. 




III. REDISTRIBUTION OF AN INPUT 
(11) ASYMMETRY 


Ref. [84| computes these terms by summing over all 
species in Fig. Both kinds of diagrams contribute 
to the asymmetry potential, while only bubble diagrams 
contribute to the thermal potential. We write the answer 
in terms of the leptons’ asymmetries, and the densities 
of the strong fluid’s conserved quantities: 


Weak processes couple leptonic and hadronic degrees 
of freedom in the primordial plasma. In this section, we 
study this coupling’s effect on lepton asymmetrie^ We 
define relevant suceptibilities in Sec. |III A[ and compute 
them over a range of temperatures in Sec. IIIB[ 




V2Gf 


' E ('*"0 - 1 


+ 2sin^ Qw)^np- 


/36{e,/r,T} 


A. Definitions and parametrization 

Let us consider the primordial plasma at temperatures 
above the quark-hadron transition temperature, Tqcd- 
The following reactions couple leptons of different flavors, 
and the quark and lepton sectors: 


-I- (l — 2 sin^ 9w) Anq 


= - 


8v^Gp 


M| 




E,, 


(12a) 

(12b) 


In the above equations, 0w is the weak mixing angle, 
and Mz/w are the masses of the weak gauge bosons. 
The symbol is a Kronecker delta, the quantities Pa 
and are net energy densities of charged and neutral 
leptons, respectively, and Aub and Anq are densities of 
the baryon number B, and electric charge Q, respectively. 
The standard model baryon number asymmetry is small 
compared to the lepton asymmetry of interest [3], hence 
can be set to zero for the purposes of this calculation. 


According to the assumptions in the first part of this 
section, the plasma starts out with a net lepton asymme¬ 
try in the mu flavor. As we showed in Sec.|lj this asymme¬ 
try is redistributed between muons and muon neutrinos. 
Moreover leptons of other flavors acquire asymmetries 
that respect Eq. Q, and the strong fluid acquires a net 
electric charge density Anq to maintain overall neutral¬ 
ity. Equation (12a) shows how the asymmetry poten¬ 
tial depends on the redistributed asymmetries, which we 
study in the ensuing section. 


Va+ P vp + a , (13a) 

Va + ot^ ^ a -I- 6, (13b) 


where a and b are quarks with charges of -1-2/3 and —1/3 
respectively. Free quarks no longer exi st at temperatures 
below TqcD, and the reactions in Eq. ( |13b ) transition to 
ones involving mesons, like Eq. ([^. 

In principle, we could study the effect of all these 
reactions on input asymmetries, but it is a daunting 
task; one that is further complicated by the quark- 
hadron transition. The following consideration of the 
relevant timescales suggests a solution. In the radi¬ 
ation dominated era, the Hubble rate \s H « 2 x 
10® s“^ g*^/^(T/GeV)^. At temperatures above the 
quark-hadron transition, the rates of reactions in (13) 
are r(r) ~ G|r® « 2 X IQi'* s"! (T/GeV)®, while the 
relevant rates below the transition are those of pion de¬ 
cays. The most important channel for the latter is the 


^ During this preparation of this manuscript, we became aware of 
Ref. 1991 . which points out the relevance of this effect to ster¬ 
ile neutrino production, and estimates it under the simplifying 
Stefan-Boltzmann approximation for free quarks. 
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muonic decay, •k'^ —>■ which is faster than the 

Hubble rate = 3.8 x 10^ s“^). Thus, a sig¬ 

nificant number of the reactions in Eq. ( |13[ ) are faster 
than the Hubble rat^Hl 

This has two primary consequences. Firstly, high re¬ 
action rates enforce kinetic equilibrium, i.e. all active 
species’ PSDs approach the Fermi-Dirac or Bose-Einstein 
forms. Secondly, forward and backward reactions are in 
chemical equilibrium, one effect of which is to equate the 
chemical potentials for both sides (the Saha equation). 
However, it has another implication - the plasma’s com¬ 
plicated internal dynamics can be abstracted into a few 
parameters or susceptibilities that completely specify its 
response to small external ‘forces’, or in this case, input 
asymmetries. All that remains is to compute the suscep¬ 
tibilities relevant to our problem. 

We now define a few useful quantities and notation. 
Given any conserved quantity F, the symbol denotes 
its chemical potential. The asymmetry Aha, in a particle 
A, is a function of its chemical potential jiA = (J-a/T. The 
quantities An^ and (la are small, and in the linearized 
limit, related by 


AfiA = XaAa, (14) 

where xa = XaIT'^ is the number-density susceptibility. 
The lepton asymmetries in the three flavors are 

ta = Aha- + Ahu^ 

= Xa-h-a-aG{e,n,T}. (15) 


The strong fluid is described by the densities of its con¬ 
served quantities: the charge and baryon-number densi¬ 
ties Atzq and Aub, respecti ve!#] The chemical equilib¬ 
rium of the reactions in Eq. ( |13| ^ implies 

- ha-- i^Q = 0, aG{e,fL,T}. (16) 


Here (lq is the chemical potential for adding a unit of 
electric charge. The conserved quantities’ densities are 
related to their chemical potentials by their susceptibili¬ 
ties: 


fAhQ\ 

{AhBj 


~QB 

A2 All 
^BQ -B 

Xll X2 



(17) 


^ The electronic channel for the pion decay, tt^ Ue is 

helicity-suppressed = 4.7 x 10^ s~^) s,nd of the 

order of the Hubble rate at temperatures T ~ 50 MeV, hence 
one might worry that leptons with electronic flavor depart from 
equilibrium. This is resolved by the observation that they are 
coupled to muonic species by other non-helicity suppressed, and 
consequently faster, reactions such as e"*" + i/g M"*” + 

O e+ + 

^ We do not follow the strangeness, S. since it is not con¬ 
served in weak reactions. Above the transitions, it disap¬ 
pears at the Cabbibo—suppressed rate Tq ~ \Vus \^~ 
10^^s~^(T/GeV)^, while below the transition the relevant rate 
is the Kaon inverse lifetime, = 8.1 x 10^ s~^. 


Equation 0, along with net charge and baryon number 
conservation, yields the constraint equations 

Aub = Xn^Aq + X2 Ab « 0, (18) 

0 = X?Aq + X?i^A b - E 

a£{e,(L,r} 

= x?Aq + x?i®A b - E • h9) 

aG{e,(L,T} 


Equations p5|, (16), (18) and are eight linear equa¬ 


tions for eight unknowns. The resulting asymmetries 
(obtained via their chemical potentials) are the ‘redis¬ 
tributed’ input lepton asymmetries Ca- 
We symbolically represent the solutions as 


hA = 


d(iA 


E - 

dC 

Q;G{e,/i,T} 


( 20 ) 


where the coefficients {d(iA/d£a) depend on the suscep¬ 
tibilities of both the leptons and the strong fluid. We 
also express the redistributed asymmetries as 


Aua = ^ 


Q;G{e,/i,T} 


OAua 

dta 




Q!G{e,/.t,r} 


, d(iA/i 

uGql 


( 21 ) 


At the temperatures of interest, the lepton susceptibili¬ 
ties are essentially given by the free particle, or Stefan- 
Boltzmann, formula: 

XAihiA) = - ^ ^ dpfh'F + ■ ( 22 ) 

In this equation, qa and toa = ruA/T are the 
spin degeneracy and mass respectively, h'p{x) = 
{d/dx){l/[exp{x) + 1]} is the derivative of the Fermi- 
Dirac distribution. The strongly interacting fluid’s sus¬ 
ceptibilities are considerably more complicated, espe¬ 
cially near the quark-hadron transition. We evalu¬ 
ate them using a number of techniques: perturbative 
quantum chromodynamics (QCD) at high temperatures, 
matching to lattice QCD results near the transition, and 
a hadron resonance gas (HRG) approximation at low 
temperatures. 


B. Susceptibilities of the strongly-interacting 
plasma 


In this section, we compute the strongly-interacting 
plasma’s susceptibilities to baryon number and electric 
charge fluctuations. The susceptibilities are the following 
derivatives of the QCD pressure 


X 2 


d'^PQCB 

dfij. 


iLx=0 


(23) 


All 


^^Aqcd 

d(ixd(iY 


fi.X,ilY=0 


(24) 
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where X,Y£ {B, Q}, fix = ^J-x/T is the chemical po¬ 
tential of the conserved charge X, and the pressure pqcd 
is given by the logarithm of the QCD partition function 
Zqct)- 


Pqcd = (25) 


where V is the volume. In Eq. (24), the off-diagonal term 
encodes the correlation between the fluctuations of 
conserved charges X and Y. Note that the susceptibili¬ 


ties in Eqs. (231 and (24) are dimensionless. 


The sterile neutrino production calculation carried in 
the present work requires knowledge of these suscepti¬ 
bilities over a broad range of temperatures, both above 
and below the quark-hadron transition. At very high 
temperatures T ^ Tqcdi the QCD pressure can be com¬ 
puted using a standard perturbative approach, while at 
intermediate temperatures T ~ Tqcd, perturbative tech¬ 
niques become inadequate; we must rely on lattice cal¬ 
culations (see e.g. | 100 lllOlj l to compute susceptibilities 
through the quark-hadron transition. At low temper¬ 
atures T < Tqcd, we compute the QCD pressure using 
the hadron resonance gas (HRG) model |1021ll03j . which 
approximates the QCD partition function as a sum over 
all known hadronic resonances. Our strategy to compute 
the susceptibility tensor over the whole required range 
of temperatures is as follows: we first separately calcu¬ 
late it both above and below the quark-hadron transition 
using either perturbative or HRG techniques, and then 
smoothly join the results with those from lattice QCD 
computations in the regions of overlap. 


1. High-Temperature Limit: Perturbative Approach 

We follow the approach of Ref. |104j to perturbatively 
compute the QCD pressure and its derivative up to order 
0{g^), where gs is the standard QCD gauge coupling 
constant. The starting point is to write the QCD pressure 
as 

PQCD = am + , ( 26 ) 

where g^ is the effective gauge coupling 

93 = 9s + ’ ( 27 ) 

and the functions are given by 
_ 2 

am = ^ + dCA , ( 28 ) 


am = -dA X! ^ 0^2 {mi,fii) [l -k 6 T 2 {rhi, fn)] 
i=i I ° 

(31»£-+2)AK.fe) 


-2miFi {ml,fLi) 


dA^A 
144 ’ 


(29) 


MS _ 22C'a 
aE7 - — 



1 


2 

3 


Ni 

E 

i=l 


2 In 


JL 

rrii 


+ F3 {ml.fii) 


(30) 


Here, dx = — 1 and Ca = stand for the gauge- 

group constants for the adjoint and fundamental repre¬ 
sentation of SU{Nc), respectively. In this work, we adopt 
the standard value of W = 3. In the above, Nf is the 
number of quark flavors, p, is the energy scale at which 
the masses and the coupling constant are evaluated (not 
to be confused with the chemical potentials), and 7 e is 
the Euler-Mascheroni constant. The functions Fi,... ,¥4 
are given in appendix [Aj 

We also need a prescription for the running of the cou¬ 
pling constant gg and of the quark masses with the energy 
scale p. As in Ref. [104] . we adopt a simple 1-loop run¬ 
ning which yields 


5s (/^) = 


2477^ 

(llGA-4TF)ln(/2/Aj^)’ 


(31) 


9Cf 

llCx-4Tp 

, (32) 

where Tf = N{/ 2 , Cp = {N^ — l)/(2fVc), and Aj^ is the 
MS renormalization scale. Here, we take pj-^f = 2 GeV. 
We follow Ref. |105j and use the criterion of minimal 
sensitivity to set the scale p 


mPp) = 


In (^ref/-^Ms) 
In (p/A^) 


-IVc+4Wf ln4 

p = Yr:Te~'^^e 22iVc-4iVf ^ 


(33) 


To compute the susceptibilities, we substitute Eq. (26) 


into Eqs. (23) and (24), remembering that the relation 


between the quark chemical potentials and those of the 
conserved baryon number and electric charges is 



(34) 


where pu and pd are the chemical potentials for up- and 
down-type quarks, respectively. We numerically evalu¬ 
ate the integrals in the functions Fi,..., F 4 . The tem¬ 
peratures relevant to sterile neutrino production are well 
below the top quark mass; we therefore adopt N{ = 5. 
We use quark masses evaluated at the reference scale p^ef 
from Ref. [T]. The only free parameter in the perturba¬ 
tive calculation is the MS renormalization scale, which 
we set so as to match with the lattice calculation in the 
regime where both approaches are valid. This entails us 
setting Aj^ = 105 MeV; the results are only logarithmi¬ 
cally dependent on the particular value we choose. 
































2. Intermediate Temperatures: Lattice Calculations 

Susceptibilities at temperatures close to the quark- 
hadron transition have been previously studied in the 
context of heavy ion collision experiments pros], where 
the crossover is signaled by fluctuations in the plasma’s 
quantum numbers |1071 llOSj . At zero chemical poten¬ 
tial, these are studied in lattice QCD by mapping to 
the expectation values of traces. Their auto- and cross¬ 
correlations are directly related to the diagonal and off- 
diagonal elements of the susceptibility matrix of Eq. 0 

[insj. 

We directly use the susceptibilities so computed at in¬ 
termediate temperatures. Specihcally, we use the results 
from the Wuppertal-Budapest (WB) lattice QCD collab¬ 
oration proo] and the HotQCD collaboration |101j . Even 
though the groups use different staggered fermion actions 
on the lattice, their results are broadly consistent with 
one another. They report the susceptibilities y®, 
and together with their estimated errors, in (2-1-1)- 
flavor QCD extrapolated to the continuum limil|j The 
lattice QCD results are in good agreement with the per¬ 
turbative calculations described above for the tempera¬ 
ture range 250 MeV < 300 MeV, above which they 
underpredict the primeval plasma’s susceptibility owing 
to the charm quark’s influence |110j . Therefore, we do 
not consider the lattice QCD calculations at tempera¬ 
tures above T > 300 MeV to avoid biasing our results. 


3. Low-Temperature Limit: Hadron Resonance Gas 


At temperatures below the quark-hadron transition, 
we model the strongly-interacting sector as a gas of 
hadronic resonances. In this HRG model, the pressure 
entering Eq. (251 is given by 


Phrg = 


vr3 


(35) 


where 


InZf = - 
InZ? = 


27r2 
VT^ 


2 Tr^ 


Y, lnZf(V,T,^Q) 

\ iG mesons 

+ 0 (E,T,/XQ,/rB) j , 

jG baryons / 

; foo __ 

-dj #p2ln(l-z,e-VP'+™'),(36) 
Jo 

roo -- 

dj dpflnil + z,e-y'p"+^^, (37) 
Jo 


where dt denotes the degeneracy factor of species f, and 
Zi is the fugacity 


— gBiAB+QiAq 


(38) 


^ The WB collaboration does not directly report X?i^, t>ut we infer 
it from their results via a change-of-basis operation. 


Here B^ and Qi are the baryon number and electric charge 
of species f, respectively. We construct the partition func¬ 
tion given in Eq. (351 by summing over all hadron res¬ 
onances with mass below 2 GeV from the particle data 
group p. We then compute the susceptibilities using 
Eqs. (23) and (24). We numerically perform the inte¬ 
grals in Eqs. (36) and (37). 

We find that the HRG results are in good agreement 
with the lattice QCD calculations for temperatures 125 
T < 150 MeV, and we smoothly match the HRG- 
derived susceptibilities to those from the lattice technique 
in this regime. 


4- Susceptibilities at all Temperatures 


We combine results from the three regimes into a sin¬ 
gle smooth susceptibility tensor, valid over the range of 
temperatures relevant to the production of sterile neu- 


trinos with masses of order 0(10 keV). Figures 4a 4b 
display the susceptibilities y®, ^^d 


4c 


and 

temperatures satisfying 10 MeV < T < 10 GeV. The 
thick solid black lines are our smooth fits to the three 
regimes, while the dashed red and cyan dotted lines are 
the HRG and perturbative results, respectively. We also 
show the results from the WB lattice QCD collaboration 
in the neighborhood of the quark-hadron transition. For 
comparison, we also show the susceptibilities computed 
in the Stefan-Boltzmann limit, i.e., assuming free quarks 
throughout and using Eq. (22). 


We observe that the HRG calculation agrees well with 
the lattice QCD result for T < 150 MeV, but systemati¬ 
cally overpredicts the susceptibilities at higher tempera¬ 
tures. The perturbative approach is consistent with the 
available lattice QCD data at T > 225 MeV, but again 
systematically overpredicts the susceptibilities at lower 
temperatures. Generally, the Stefan-Boltzmann approxi¬ 
mation overestimates the susceptibilities by a factor of or¬ 
der unity, except near the quark-hadron transition. Inter¬ 
estingly, we observe an accidental cancellation in the off- 
diagonal susceptibility, y?i^, in the (2-|-l)-flavor model 
which does not appear in the W = 5 theory. This arises 
because the sum of the electric charges of the up, down, 
and strange quarks exactly vanishes. Hence, we expect 
y^y —> 0 for temperatures above the strange quark mass 
in the (2-|-l)-flavor model. In the Vf = 5 model however, 
the charm quark becomes rapidly important at T > 300 
MeV, leading to a sharp turnover in y^y near this tem¬ 
perature. 

Given a set of infinitesimal lepton asymmetries, we 
solve for the chemical potentials using the above suscep¬ 
tibilities in Eqs. (15), (16), (18) and (fl^. We obtain the 


redistributed asymmetries in all the constituent species 
by using these chemical potentials, along with the ap¬ 
propriate susceptibilities in Eq. ( |^ . Figure [4d| plots the 
redistributed asymmetries for an infinitesimal input mu 
leptonic asymmetry. We note the following features: 
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(a) electric charge subsceptibility 





FIG. 4: Panels (a)-(c); Components of the quadratic susceptibility tensor for the primordial plasma’s electric charge and 
baryon number. In all panels, the thick solid black line shows our smooth fit used in the computation of sterile neutrino 
production. At low temperatures, we illustrate the HRG results with dashed red lines, while the high-temperature 
perturbative results are shown with dotted cyan lines. We also show the results from the WB lattice QCD collaboration [100] 
with green errorbars. For comparison, we also display the Stefan-Boltzmann approximation to the susceptibilities assuming 
free quarks at all temperatures. Panel (d): Effective populations of all leptonic degrees of freedom after the redistribution of 
an infinitesimal mu leptonic asymmetry at all temperatures. 


1. At temperatures T > 2 GeV, the redistribution is 
efficient and ~ 60% of the mu leptonic asymmetry 
ends up in the muons. All the charged leptons are 
effectively massless at this epoch, hence the popu¬ 
lations of the electron and tau flavors are identical. 


2. The rise in the mu and tau lepton populations 
above temperatures of ~ 25 MeV and 300 MeV 
reflects, in part, the rise in their particle number 
susceptibilities as the temperature becomes com¬ 
parable to their masses [see Eq. (22|]. However, 
the largest contribution to the former is from the 
disappearance of the hadronic degrees of freedom 
below the quark-hadron transition, and the associ¬ 
ated drop in the strongly interacting fluid’s suscep¬ 
tibilities. 


3. The ‘kink’ in all the redistributed asymmetries 


close to temperatures T ~ 170 MeV is a signa¬ 
ture of the sharp change in the strongly interacting 
fluid’s susceptibilities at the quark-hadron transi¬ 


tion [see Figs. 4a 4b and 4c 


4. At lower temperatures, T < 30 MeV, the redis¬ 
tribution is inefficient and most of the asymmetry 
ends up in the muon neutrinos. Moreover, the elec¬ 
tron neutrino and the muon have identical (small) 
populations. This is characteristic of inelastic neu¬ 
trino scattering, -\- e~ —>■ , which is 

the most important channel at these temperatures 
(the hadronic susceptibilities are negligible at this 
epoch). 


These redistributed asymmetries impact sterile neu¬ 
trino production via the asymmetry potential, . 
Equation (12a) expresses this potential in terms of the 
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text of neutrino decoupling, active-active neutrino os¬ 
cillations and supernova calculations d nn [m]. In 
particular, Ref. m lists a number of relevant matrix 
elements. Our calculations apply to earlier epochs, with 
a larger number of reactions due to the population of 
hadronic species above the quark-hadron transition. 

Early work on sterile neutrino production used simple 
prescriptions for the resultant increase in reaction rates 
[MES]. Recent work in Refs. diiEni provides a theoreti¬ 
cal framework to include particle masses and statistics in 
the neutrino opacity calculation, and formalism for loop 
corrections. We include a number of additional contri¬ 
butions to the neutrino opacities that are significant at 
the temperatures relevant to sterile neutrino production. 
We adopt the following simplifying assumptions: 


FIG. 5: Asymmetry potential per unit unsealed mu 
lepton asymmetry: Solid line shows the effect of 
redistribution using a combination of perturbative QCD, 
lattice calculations and the HRG approximation. Dashed 
line shows the result using Stefan-Boltzmann approximation 
for free quarks. The value is constant (= 2^/2 = 2.83) in the 
absence of redistribution [see Eq. ( |12a| )]. The redistribution 
is a ~ 10% correction at the most relevant temperatures for 
sterile neutrino production (T > 100 MeV). 


asymmetries in the populations of the individual charged 
and neutral leptons, along with those in the charge and 
baryon number of the strongly interacting fluid. As ear¬ 
lier, for an infinitesimal input mu leptonic asymmetry, 
the individual asymmetries are formally represented by 
the functions in Eq. (21); the solutions for the charged 
and neutral leptons are as plotted in Fig. We ob¬ 
tain the electric charge density of the strongly inter¬ 
acting fluid, Auq, using net electric charge neutral¬ 
ity, i.e. Eqs. (17) and (19). Tables of susceptibili¬ 
ties, along with the functions in Eq. (21) at a num¬ 
ber of temperatures from 10 GeV down to 10 MeV can 
be found at https : //github. com/ntveem/sterile-dm/ 
tree/master/data/tables, 

Figure shows the potential per unit physical p, lepton 
asymmetry using these solutions; this quantity is con¬ 
stant and equals 2\f2 = 2.83 in the absence of redistri¬ 
bution. As shown in the figure, asymmetry redistribu¬ 
tion corrects the potential at the ten-percent level above 
temperatures T > 100 MeV, which is where the bulk 
of the sterile neutrinos are produced. This contribution 
changes the resonant momenta, and the resultant sterile 
neutrino dark matter’s phase-space densities; we explore 
this further in Sec. IVl 


IV. NEUTRINO OPACITY 


1. We neglect small asymmetries in the participating 
species’ populations (as for the thermal potential). 
This is justified since the scattering rates are non¬ 
zero even in a CP symmetric plasma. Moreover, 
we assume thermal and kinetic equilibrium, due 
to which the populations of all active species are 
Fermi-Dirac/Bose-Einstein distributions. 

2. We integrate out the massive gauge bosons, Z and 

and approximate the weak interaction by a 
four-fermion contact term. Consequently, the re¬ 
actions separate into leptonic and hadronic pro¬ 
cesses, depending on the species involved. More¬ 
over, we neglect the thermal populations of and 
. These steps are valid at low temperatures and 
momentum transfers, i.e., T^sjtju <C My/±/zo ~ 
80 GeV. We operate in the temperature and energy 
ranges 

10 MeV < T < 10 GeV, (39) 

<E^JT <20. (40) 

The approximation fails at the higher energies at 
the upper end of the temperature range. However, 
as we see in Sec.|Vj the bulk of the sterile neutrinos 
are produced at lower temperatures. 

3. We assume incoming and outgoing particles to be 
non-interacting within two limits - below and above 
the quark-hadron transition (see §3.3 of Ref. |113) ). 
Below the transition, we include hadronic channels 
with pseudoscalar and vector meson^ and neglect 
the small population of baryons. Above the transi¬ 
tion, we include reactions with free quarks, i.e. we 
neglect the strong coupling constant. 

This approximation fails at temperatures T ~ 
Tqcd IRMj. We show opacities interpolated 
through the transition using a few prescriptions, 
whose consequences for sterile neutrino production 
we explore in Sec. [V] 


In this section, we outline our calculation of muon neu¬ 
trino opacities in the early universe. Initial work in this 
area focused on reactions involving leptons, in the con- 


We also include quark production in s-channel reactions at high 
CM energies. See Sec. |B I'^ for details. 
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TABLE I: Two-particle to two-particle leptonic reactions 
contributing to the muon neutrino opacity: antineutrinos are 
similar. Symbol a runs over the other leptonic flavors, i.e. 

a € {e,r}, and stands for u/i/. 


s-channel t-channel mixed 


- 

->• I'a + Q:'*’ 

^IJ, -i- - 


+ P M ■ 

Vfi + V ^ 

-1- - 

->■ I'a a 

-i- V a - 

+ V a 

P^ -1- - 

->■ 

-1- - 

> 0(~ + cx'^ 

v^+a~- 

+ Va - 

+ Va 

> + «■•■ 

P^ -1- Pm - 

> 


The collision integral for a massless muon neutrino is 


+ (41) 

where T and are the interaction rate (opacity) and 
PSD, respectively (This expression satisfies detailed bal¬ 
ance, see assumption #1). The interaction rate is given 
by a sum over all reactions that consume the muon neu¬ 
trino [see Eq. 0]. 

It is useful to define the scaled interaction rate 




GlT^E^ 


(42) 


In the limit where all the particles involved are relativis¬ 
tic, weak cross-sections are proportional to the squared 
energy in the CM reference frame. If we ignore parti¬ 
cle statistics, reaction rates follow the scaling of Eq. 0, 
hence the scaled rate is proportional to the number of 
relativistic degrees of freedom involved [M]. We present 
the scaled rates in the rest of this section, in order to 
contrast our results with this intuition. 

In the rest of the section, we enumerate reactions con¬ 
tributing to the opacity and present final rates under the 
above approximations. We elaborate the calculation of 
matrix elements in Appendix [B| 

Table |l] lists the purely leptonic two-particle to two- 
particle reactions contributing to the muon neutrino 
opacity. It is similar to the list in Ref. albeit 

with tau leptons, which are important at temperatures 
T > 400 MeV. We also include ‘three-body fusions’, since 
they arise at the same level of approximation. These are 
generated by omitting in turn the products in the reac¬ 
tions of Table H) adding their charge conjugates to the re¬ 
actants, subject to the constraint that the product’s rest 
mass is strictly greater than the sum of the reactants’. An 
example is tau lepton production via -|- /i+ -|- 11 "+. 

Figure shows the leptonic contribution to the muon 
neutrino opacity at a temperature T = 100 MeV, using 
the matrix elements for reactions in Table [ij and related 
three-body fusions. For convenience, we only show re¬ 
actions in the top hve at any particular momentum bin. 
In the numerical implementation, we evaluate the dimen¬ 
sionless quantity /G'^T^ (proportional to the un¬ 

sealed rates) to an accuracy of 10“® after simplifying the 
collision integrals in Eq. (B7) and (B25). 



U^Uee~ fJ.~ 

{-) (-) 
U l/fj, u 

U^T~ —>• fl~UT 

VftOe 


FIG. 6: Scaled muon neutrino opacities through leptonic 
reactions, vs energy at T = 100 MeV. Only reactions in the 
top five at any particular momentum bin are shown. The 

symbol stands for n/u. 


The quark-hadron transition considerably complicates 
the hadronic reactions. We appeal to assumption # 3 
and evaluate their rates in two limits: at low and high 
temperatures, i.e. T < Tqcd and T > Tqcd respectively. 
Table [n] lists the hadronic two-particle to two-particle re¬ 
actions contributing to the muon neutrino opacity. At 
high temperatures, we list all reactions involving free 
quarks, while at temperatures T < Tqcd , we assume that 
all incoming hadronic degrees of freedom belong to the 
pseudoscalar meson nonet (the heaviest member of which 
is the eta meson, with a mass of rurj = 547.8 MeV). As 
in the leptonic case, we also include three-body fusions 
involving pseudoscalar mesons or quarks. Examples are 
and charm quark production via -\- ■k~ —>■ 

and s ^ c, respectively. 


A complication is that at low temperatures, free quark 
and parton currents contribute to the initial and final 
states for large momentum transfer in the t- and s- 
channel respectively (see Appendix B 1 a). For s-channel 
reactions, we thus treat individual meson resonances for 
center of mass energies < 1 GeV, and use the free quark 
model for inclusive cross sections at > 1 GeV. Also im¬ 
portant are ‘two-body fusions’, i.e. reactions with two 
particles in the initial state and one in the final state, 
with the latter being a pseudoscalar or vector meson. 
Table |V| in Appendix |B 2| lists all such reactions included 
in our opacities. 


Figure shows the hadronic contribution to the muon 
neutrino opacity at low and high temperatures, using the 
matrix elements for two-particle to two-particle reactions 
in Table |0 the associated three-body fusions, and two- 
body fusions in Table |0 of Appendix |B 2[ As earlier, we 
only show reactions in the top five at any momentum bin; 
the numerical implementation of the first two classes is 
unchanged. 
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TABLE II: Two-particle to two-particle hadronic reactions 
contributing to the /r neutrino opacity: antineutrinos are 
similar. Symbols a and b run over quarks with charge 

-|-(2/3)e, and —(l/3)e, respectively, and stands for a/a. 


s-channel 


t-channel 


T > Tqcd 


-1- /r+ —>• a -1- 6 


1 V 1 

“h CL a 

1/ fj, u ^ ^ a d 


i^M + ^ —>■ I'm + ^ 

Vfx + Ufx — y b b 


Ufj, b ^ fj, a 



d ^ fjj -\-b 


T < Tqce0 

I'm + 7r=‘= —1 I'm -f tt* 

^ 7r° 


Uf, + p+ ^ K+ +K0 


I'm + ^ I'm + 

Ufj, + —>■ tt’*' + 


I'm 

—>■ + TT*^ 


n^ + K- ^ p- +K^ 

^ + V 


I'm + 7“” 

Vfj, ^ TV'^ tt” 


Vpi + K —>■ fjj + 7V^ 

+ I'm 


I'm + K~ p~ -\-p 

I'm + 7r° —>■ -1- 7r+ 

+ ^ p- + K+ 

K° ^ p~ ->r 7r+ 
I'm -t- 7T^ —y p 4- 

“'m + “7 —^ + K* 


Input neutrinos can produce quarks (s-channel) or probe 
mesons’ quark content (t-channel) at T < Tcic t) for large 
momentum transfer in Z^, . See Appendix |B 1 alfor details. 


Figure shows the total opacities with muon neutrino 
energy at temperatures of 100 MeV and 2 GeV. We note 
a few salient features of these rates. 

Firstly, we note that the leptonic and hadronic two- 
particle to two-particle reaction rates approach the scal¬ 
ing of Eq. (§ at large energies; the downturn at lower 
energies is due to Pauli blocking. 

Secondly, both sets of rates exhibit a rise at low en¬ 
ergies, which reflects non-zero limiting values of the un¬ 
sealed rates. This is due to the behavior of cross sections 
for inelastic collisions involving massive particles, such 
as the three body collision e~ Pg or the 

scattering process -P- Vg-\-e'^. We illustrate this 

by calculating the cross-section for the latter, while ne¬ 
glecting the positron’s rest mass and Pauli blocking for 
simplicity. The squared and spin-summed/averaged ma¬ 
trix element for this process is 

(|7W|2) = 128G| ■pg+) {p^+ •Pj.J . (43) 

In the limit of zero neutrino energy 

{p^L+ ■Pu^)= -ml/2 + O , (44) 

which implies that modulus squared is 

(|M|2) = -MGlml ■Pg+)+0 (G.J . (45) 


- fi 7 r+ 

- - i/^p+7r“ -> 

. > /A“ 7 r^ 

- ■ - K + 

- 7 r+ 

- 

. Uii7r~ —f 

- • - 

- - -> 7 r+ 7 r° 

■ - 

_ p® 

. ud 

- I'fii'fi —> dd 


10 -“ 10 -=* 10-^2 10-1 10 " 10 “ 
= E^JT 

(a) T = 100 MeV 



|3 


— i/^c 

■. - Uf^n+s c 

. Uf^c fj,~d 

- u^H^d c 

- - Ufj^d —)• fj,~u 

ud 

UfiU —V fj.~d 
- Uf^S ii~c 


10 -“ 10 -" 10 -" 10 -“ 10 " 10 “ 
E^^ = E^JT 

(b) T = 2 GeV 



FIG. 7: Scaled muon neutrino opacities through hadronic 
reactions, vs energy: panels (a) and (b) show rates at 
T = 100 MeV and 2 GeV, respectively. Only reactions in the 
top five at any particular momentum bin are shown. 


Hence the cross section for the p—neutrino, integrated 
over outgoing particles’ directions, is 

Ciml 

a., = ^:^+0{E,^). (46) 

Such non-zero limiting values are responsible for the rise 
in the scaled rates for soft neutrinos in Figs. |^and[^ 
Thirdly, the hadronic opacities at low temperatures, 
i.e. T < Tqcd, exhibit a series of peaks. These are sig¬ 
natures of two-body fusions, which are broad resonances 
in the propagators of the weak gauge bosons. These in¬ 
clude the production of pseudoscalar mesons (e.g. pion 
production via -|- 7r+) and vector mesons (e.g. 

production via Pfj, ^ p°). In the total opacities 
of Fig. 1^ the former is visible as a peak at intermedi¬ 
ate momenta, while the latter are smeared out at large 
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T [MeV] 


(a) rates vs energy 


(b) rates vs temperature 



login (K = (-E.n = E,JT) 


(c) Tc = 250 MeV 


(d) Tc = 1000 MeV 


FIG. 8: Scaled muon neutrino opacities for a range of energies and temperatures: panel (a) shows total, leptonic, and 
hadronic opacities vs energy at T = 100 MeV and 2 GeV. Panel (b) shows opacities at high energies {Ei^^ — 20T) vs 
temperature: black lines are two interpolations through Tqcd- They are cubic splines labeled by their cutoff temperature, T^, 
as defined in the text. Golored lines are numbers of relativistic degrees of freedom: g*, KT under assumption ^ 3, i.e., that of 
Kolb and Turner m with Tqcd = 210 MeV, chosen to match Ref. mM whose results are g*,LS(2006). Panels (c) and (d) 
show interpolated opacities vs energy and temperature. Blue dashed lines mark ranges where the two values of g* differ by 
more than 10%, and are a rough guide to where these rates can be trusted. Red dotted lines mark the cutoff temperature. 


momenta. 


Finally, we observe from Fig. 8a that the total opacities 
at high energies exhibit a jump as the temperature passes 
through Tqcd ■ This is due to the increase in the number 
of hadronic degrees of freedom, as evidenced by the sizes 
of the jumps in hadronic- and leptonic contributions (the 
latter due to the tau lepton turning on). 


This is shown clearly in Fig.|8b[ which shows the scaled 
muon neutrino opacities at high energies for a range of 
temperatures. Note that these rates assume that the 


hadronic species above and below the transition are free 
quarks and mesons respectively (assumption 3 in our 
list above). For comparison, the figure shows the num¬ 
ber of relativistic degrees of freedom, g,, both under this 
assumption and from Ref. [104] . which implements the 
running of the strong coupling constant. We note the 
significant deviation close to the quark-hadron transition 
(Tqcd = 210 MeV in the lattice calculations underlying 
Ref. [TTH] !. 

Motivated by this, we explore two methods of interpo- 
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lating opacities through the quark-hadron transition. In 
each of them, we choose a cutoff temperature, Tc, above 
which we use the free quark results, and use a cubic spline 
interpolation in between. We emphasize that this is not 
physically motivated; the actual rates, and their matrix 
elements, need to incorporate the strong coupling con¬ 
stant and its running. The figure shows interpolations 
with Tc = 250 MeV and 1000 MeV, which we expect to 
bracket the range of rates. 

With this caveat, Figs. |8c| and [8d| show interpolated /i 
neutrino opacities for a range of energies and tempera¬ 
tures, with Tc = 250 MeV and 1000 MeV, respectively. 
In the rest of the paper, we use these scattering rates 
and the potentials defined in Sec. to study sterile neu¬ 
trino production via oscillation in the early universe. We 
use both interpolations through the quark-hadron tran¬ 
sition in order to illustrate the results’ sensitivity to the 
scattering rates. 


V. STERILE NEUTRINO PRODUCTION 


In this section, we incorporate the standard model cal¬ 
culations of Sections m and IIVI into the sterile neutrino 
production mechanism, whose broad outline we provided 

in Sec. [TTl 

We evolve the sterile neutrino and anti-neutrino PSDs, 
/I'sip) /pa(p)) using the Boltzmann equation of 

Eq. (|^. We use the primordial plasma’s temperature 
T as a, clock, and numerically integrate a thousand log¬ 
arithmically spaced Lagrangian momentum bins from a 
temperature of 10 GeV down to 10 MeV. For the mod¬ 
els illustrated by stars in Fig. [T| the vast majority of the 
production happens between these temperatures. We use 
the muon neutrino opacities derived in Sec. |1V| and pro¬ 
vide results using the two interpolation scheme presented 
in Fig. which bracket the range of uncertai nties due 
to the quark-hadron transition. We use Eq. (12b) for 
the thermal potential and the results presented in 
Fig. ^for the asymmetry potential incorporating the 
redi^ibution of Sec. nn 

In order to close the system of equations, we also need 
the evolution of the plasma temperature T and mu lep- 
tonic asymmetry with coordinate time t. Before dis¬ 
cussing the details of the sterile neutrino production, we 
briefly review these two relations. 


A. Time-temperature relation 

In this subsection, we derive the time-temperature re¬ 
lationship prior to the epoch of weak decoupling. The 
Hubble rate, H, is 


where a is the scale factor, mp = 1.2 x 10^® GeV is the 
Planck mass, and psM and p,j^ are energy densities in 
standard model particles and sterile neutrinos, respec¬ 
tively. The latter is given by an integral over the PSDs, 
Pi'. = (l/27r^)/p^dpA/p2 + f^Ap)]- During 

Hubble expansion from a to a+Sa: a) the sterile neutrino 
PSDs evolve to + 5f^^/c^{p) due to a combina¬ 

tion of mixing with muon neutrinos, and their momentum 
redshifting as pa = const, b) due to large neutrino opaci¬ 
ties, all active species maintain equilibrium PSDs with a 
common temperature, whose evolution is affected by the 
production of sterile neutrinos. 

The continuity equation for the total stress-energy ten¬ 
sor is 

~ ~ dr ^ ’ 

(48) 


where Psu/v. are SM and sterile neutrino pressures, re¬ 
spectively. The sterile energy density evolves according 
to 


dpu, dpv 


dT 5 In a dT 


dlna ^ 


dt 


dt dT' 


(49) 


The two terms on the right-hand side are the free- 
streaming and oscillation contributions, respectively. 


dpu, 

d\na 

dPi'. 

dt 


— ~^iPi'. + 
f dp p^ 


(50) 


271^ 


-Vp' 


^ (p) -h (p)] . (51) 


We substitute Eqs. ( |49[ ) and ( |50| ) into Eq. ( |48[ ) and solve 
for the relation between the scale factor and temperature 


dlno 


/ dpsM 

I dr 


dpv^ dt 

DTdr 


(dSM + ^’sm) ^ • (52) 


Substituting Eq. (471, we obtain the time-temperature 
relatioi0 


dT _ 3jd[psM + Psm] + {dpvjdt) 
dt dpsu/dT 


(53) 


Defining the number of SM relativistic degrees of freedom 
for the energy and entropy densities via 


SSM = 


Psm 


Psm = 

^ (55) 


we have the final form of the time-temperature relation 
dT AHg^^sT^ + {30/TT^){dp^Jdt) 


dt 


d[g,T^]/dT 


(56) 


We use numbers of relativistic degrees of freedom p* and 
p» s from Ref. |104j in our numerical implementation. 


d 

dt 


In a = id = 



(PSM + Pi-J) 


(47) 


® We note that we correct here an error introduced in Ref. m- 
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FIG. 9: We illustrate the temperature-evolution of the 
sterile neutrino’s PSD for the central model of Fig. with 
(msjSin^ 26) = (7.1 keV, 4 x 10“^^). Solid and dashed lines 
distinguish results with neutrino opacities from Fig. |8c|and 
8d respectively. 


B. Time-evolution of asymmetry 


The temperature-scaled muon asymmetry, evolves 
both from the depletion of relativistic degrees of freedom 
due to annihilations and from the production of sterile 
neutrinos. There are subtleties in dealing with the lat¬ 
ter in the case of resonant production |114| . but for the 
semi-classical approach outlined in Sec. 


IB[ we can write 


down the contribution in terms of the evolution of sterile 
neutrino PSD. Keeping in mind the definition of the lep¬ 
ton asymmetry in Eq. Q, the asymmetry evolution due 
to both contributions together is 


dCf, 

dt 


- ^ 
dt , 

= -3 


dpp^ 


27r^ 


(P) - /p;. (P) + 2/^- (p) - 2/^+ (p)] 


H 


dlnT 

dt 




dpp^ d 
277^ dt 


[UM - UM] : 

(57) 


where the symbol p is the temperature-scaled momen¬ 
tum, p = p/T. The hrst term in the square bracket 
in the last line above can be evaluated with the help of 
Eqs. (47) and ( |5^, while the second term can be eval¬ 
uated using Eq. (1^. Our large number of momentum 
bins (1000) allows us to use spline integration at every 
time step in order to perform the momentum integrals in 
Eqs. (57) and (51). We set up our Lagrangian momen¬ 
tum bins such that 5 x 10“^ < p/T < 20 at temperature 
T = 10 GeV. We have checked that this range is more 
than sufficient to accurately capture the most relevant 
range of the sterile neutrino PSDs. 


C. Resonant Production 

As described in Sec. |llj the presence of a lepton asym¬ 
metry leads to a resonant production of sterile neutrinos 
with specific momenta. Through Eq. ([^, the resonant 
momenta at a particular temperature satisfy 


A(p) cos 26-V^- W^(p) = 0. 


(58) 


Substituting the definition of A(p) and the potentials 
from Eq. (12), we obtain 


ml 

2p 


dV^ „ dV^'-{p) _ ^ 

dp 


(59) 


There are two roots, i.e. two momenta resonant at any 
temperature |18j . Consideration of the terms’ approx¬ 
imate temperature scaling shows that each scaled root 
(p = p/T) sweeps to larger values at lower temperatures 
(ignoring changes in the numbers of relativistic degrees 
of freedom). This is reflected in Fig. [9 which shows 
the sterile neutrino PSD’s evolution with temperature 
for the central model in Fig. [^with m-s = 7.1 keV and 
sin^ 20 = 4 X 10“^^. We observe that most of the neutri¬ 
nos are produced at the lower resonance and at tempera- 
ture s close to Tqcd- This is illustrated by Figs. 1 10a] and 
10c which show the evolution of the entropy-sca. tecQp 
lepton asymmetry and the net sterile neutrino and an¬ 
tineutrino density for the range of models marked by 
stars in Fig. The latter is also sensitive to thermal 
(nonresonant) production, which operates at all temper¬ 
atures, but is subdominant for the mixing angles of in¬ 
terest. 


Figures. 10b and lOd show the sterile neutrino and 
antineutrino PSDs at T = 10 MeV for these models. 
We note that the sterile antineutrinos are produced off- 
resonance for the positive lepton asymmetries we con¬ 
sider here, and their abundance is thus significantly sup¬ 
pressed compared to that of the sterile neutrinos. Solid 
and dashed lines in Fig.[^show results for the two inter¬ 
polations of the /i neutrino opacities through Tqcd pre¬ 
sented in Fig. 1^ which differ in the temperature range 
150 MeV < T < 1 GeV. For small values of the mixing 
angle, we observe that there is little difference between 
the PSDs computed using our two different interpola¬ 
tion schemes for the neutrino opacity. For these models, 
most of the production happens at temperature below the 
quark-hadron transition where our two opacity approx¬ 
imation schemes are essentially the same, hence leading 
to similar PSDs. As the mixing angle is increased, the 
pro ducti on is pushed toward higher temperatures (see 
Fig. 10c) where the difference between our two interpo¬ 


lation schemes is greater, leading to a larger uncertainties 
in the final PSDs. 


^ We show this scaling rather than the one with temperature, since 
it is conserved through epochs of annihilation. 
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(a) mu lepton asymmetry (b) sterile neutrino PSDs at T = 10 MeV 



(c) net sterile density (d) sterile antineutrino PSDs at T = 10 MeV 


FIG. 10: Sterile neutrino production mechanism: Panels (a) and (c) show the entropy scaled mu lepton asymmetry and the 
net sterile number density with temperature. For each model with a given mass and mixing angle, the mu lepton asymmetry 
at high temperatures is tuned by hand to produce the right relic abundance. Panels (b) and (d) show sterile neutrino and 
antineutrino PSDs, respectively, at T = 10 MeV. Colors differentiate models in Fig. and solid and dashed lines distinguish 
results with neutrino opacities from Fig. |8c| and |8d| respectively. Note the different numerical factors multiplying the y-axis of 
panels (b) and (d). The dotted line in panel (b) is a massless Fermi-Dirac distribution with degeneracy g = 0.003. 


Table |TT^ lists parameters describing the production 
and final sterile neutrino DM PSDs for the models 
marked in Fig.[^ Also provided are the ranges for differ¬ 
ent interpolated /i neutrino opacities through the quark- 
hadron transition as in Fig.[^ Note that the sterile PSDs 
in Figs. |10b| and |10d] are non-thermal; we show the mean 
momentum (p/T) relative to the active neutrino temper¬ 
ature scale. 

A key element to take away from Table |III| and 
Figs. |10b| and |10d] is that the ‘warmer’ models with larger 
values of (p/T) are less sensitive to the uncertainty in the 
quark-hadron transition. This is important since these 
warmer models can be most easily constrained by small- 
scale structure formation. Therefore, uncertainties in the 
strong plasma near Tqcd are unlikely to affect the ro¬ 
bustness of the these constraints. 


VI. TRANSFER FUNCTIONS FOR MATTER 
FLUCTUATIONS 

In this section, we study the effect of sterile neutri¬ 
nos on the growth of density fluctuations in the early 
universe. We focus on the lepton asymmetry-driven 
mechanism outlined in Sec. |Tlj and on modes of the 
matter distribution with co-moving wavenumbers k € 
[1,100] /iMpc”^. These scales are probed by the Lyman- 
a forest in quasar spectra (see m and references 
therein), and populations of dwarf galaxies in the Local 
Group (see [11611117j and references therein). All these 
scales enter the horizon after the redshift 2 ;h — 4 x 10^, 
when the temperature of the photon-baryon plasma is 
T ~ 10 keV. The sterile neutrino models shown in Fig. 
cease to be produced below temperatures T ^ 100 MeV; 


















































TABLE III: Parameters for the models marked in Fig. 
with rris = 7.1 keV and Sloivih^ = 0.119 [2]. The ranges 
displayed in the three last columns account for the 
uncertainties in the neutrino opacities near the 
quark-hadron transition. 


sin^ 26 

(T^/ssm); 

(T^/sSM)f 

{p/T}R 


at T — 

10 GeV 

at T = 

10 MeV 



xl0““ 

1 

O 

X 

1 

O 

X 



0.800 

13.0 

- 13.1 

6.95 

- 7.03 

2.60 - 

2.61 

1.104 

10.80 

- 10.88 

4.74 

- 4.81 

2.45 - 

2.47 

1.523 

9.57 

- 9.64 

3.51 

- 3.58 

2.28 - 

2.32 

2.101 

8.81 

- 8.88 

2.76 

- 2.83 

2.12 - 

2.16 

2.899 

8.32 

- 8.39 

2.27 

- 2.34 

1.95 - 

2.01 

4.000 

7.96 

- 8.03 

1.93 

- 2.00 

1.80 - 

1.87 

5.519 

7.69 

- 7.76 

1.68 

- 1.74 

1.66 - 

1.74 

7.615 

7.45 

- 7.53 

1.47 

- 1.54 

1.53 - 

1.62 

10.506 

7.20 

- 7.29 

1.28 

- 1.36 

1.43 - 

1.52 

14.496 

6.95 

- 7.05 

1.09 

- 1.18 

1.35 - 

1.44 

20.000 

6.7 

- 6.8 

0.9 

- 1.0 

1.29 - 

1.38 


^ The sterile DM distributions are non-thermal; we compute 
{p/T) using the active neutrino temperature. Below the epoch 
of annihilation, the latter is related to the CMB 
temperature by the factor (4/11)^/® = 0.714. We note that for 
a Fermi-Dirac distribution (p/T) ~ 3.15 

hence we can assume they are essentially collisionless in 
this section. 

The main effect of such a collisionless component on 
matter fluctuations is suppression due to free-streaming 
in the epochs where it is relativistic [mninni . Previ¬ 
ous works extensively studied this in the context of warm 
and/or neutrino DM models (see Refs. [12011121] and ref¬ 
erences therein), and identified the characteristic scales 
at which the suppression set in as a function of the neu¬ 
trinos’ mass and mean momentum [18] . 

In order to obtain the suppression’s detailed form, we 
need to incorporate the PSDs of the sterile neutrinos 
and antineutrinos into the Boltzmann equation for the 
DM component. This entails solving a perturbed form 
of Eq. (|^, with additional terms due to inhomogeneities, 
but without the source (production) terms. The scales 
of interest are non-linear in the current epoch, but we 
only provide the linear transfer functions at z = 0, which 
can be used as initial conditions for cosmological N-body 
simulations. 

We use the publicly available CLASS solver [86]to in¬ 
tegrate the perturbed linear Boltzmann equatioij^ We 
initiate the solver with the Planck background param¬ 
eters |3], except with the CDM component replaced by 


Our choice was motivated by the availability of well-documented 
modules to deal with non-cold relics. We have checked our results 
against those obtained from a modified version of the publicly 
available CAMB solver mi- 
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FIG. 11: Suppression of the transfer functions of overall 
density fluctuations relative to the ACDM ones for sterile 
neutrino models in Fig. as a function of wavenumber. 
Dashed and dotted lines show results for the interpolated fj, 
neutrino scattering rates of Figs. |8c|and |8d| respectively. 

The solid black line is the numerical fit for a thermal warm 
DM transfer function as given in Ref. mM- 

collisionless components with PSDs as shown in Figs. jlObj 
and |10d| Since we are interested in the detailed shape 
of the transfer function, we turn off the default fluid ap¬ 
proximation for non-cold relics [123j . 

Figure[^shows the resulting suppression as a function 
of the comoving wavenumber. We illustrate the sup¬ 
pression in the fluctuations’ transfer functions relative 
to their values in ACDM. Also shown is the commonly- 
used fit to a thermal warm DM transfer function given 
in Refs. [531 1154] with an ‘equivalent thermal mass’ of 
mth = 2.2 keV; fits for models marked with stars in Fig.[^ 
have a range of 1.6 to 3.2 keV. 

However, the strong difference in shape with the ther¬ 
mal WDM transfer function warrants use of the exact 
sterile neutrino dark matter transfer functions. The ther¬ 
mal warm DM PSDs relevant to the fit are rescaled ver¬ 
sions of the Fermi-Dirac distribution; as can be seen from 
Fig. |10b[ the resonantly-produced DM’s PSD has an ex¬ 
cess at low momenta that cannot be reproduced by such 
a rescaling. Hence, our DM transfer functions do not ex¬ 
hibit the fits’ steep ^ dependence at large wavenum¬ 
bers and the resultant severe suppression of power on 
small scales. This indicates that the models considered 
in the present work are more likely to be in agreement 
with small-scale structure formation constraints, as re¬ 
cently pointed out in Refs. [57M51135] . 


VII. DISCUSSION AND CONCLUSIONS 

Sterile neutrinos are a well-motivated extension of the 
standard model of particle physics, and offer a promising 
candidate for the inferred DM population of the Uni- 


(28) X 10^ 
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verse. In this paper, we performed a detailed study of 
the resonant production of sterile neutrinos with masses 
and mixing angles relevant to the recent X-ray excess. 
In doing so, we explored the rich phenomenology asso¬ 
ciated with the active neutrinos’ weak interaction with 
the primordial plasma. These interactions efficiently re¬ 
distribute primordial lepton asymmetries among all the 
available degrees of freedom, and impact the tempera¬ 
ture and momentum dependence of neutrino opacities. 
We incorporated these effects into the sterile produc¬ 
tion calculation, corrected and extended the existing nu¬ 
merical implementation, and obtained revised DM phase 
space densities. We finally computed transfer functions 
for fluctuations in the matter density, which can be used 
as starting points for 7V-body simulations of cosmological 
structure formation. 


For the parameters relevant to the X-ray excess, reso¬ 
nant sterile neutrino production coincidentally occurs in 
the vicinity of the quark-hadron transition (see Fig. 10c). 
Strongly interacting degrees of freedom affect the pro¬ 
duction in two ways: a) they influence both asymmetry 
redistribution and neutrino opacities through their inter¬ 
action with the weak gauge bosons {Z and W^), and b) 
the transition from free quarks to hadrons at Tqcd influ¬ 
ences the time-temperature relation [Eq. (56l]. We now 
consider the robustness of each of these elements to the 
remaining uncertainties in the quark-hadron transition. 

The asymmetry redistribution among the strongly- 
interacting degrees of freedom depends on the suscepti¬ 
bility of the quark-hadron plasma to baryon number and 
electric charge fluctuations. At high temperatures, we 
use tree-level perturbative QCD to compute the suscep¬ 
tibilities. There are uncertainties concerning the exact 
values of the quark masses, loop corrections, and the ex¬ 
act implementation of the MS renormalization scheme. 
We expect these to have little effect on the final ster¬ 
ile neutrino PSDs since the bulk of the production oc¬ 
curs at lower temperatures, where the lattice QCD- and 
HRG-derived susceptibilities are most relevant. Thus, 
uncertainties in the asymmetry redistribution are likely 
dominated by systematic errors in the lattice calculations 
m. measurement errors in the hadronic resonances’ 
masses, and inaccuracies inherent in the HRG approach 
near the quark-hadron transition. Our confidence in the 
fit we use in this work is bolstered by the facts that a) 
an independent lattice QCD calculation |101] find very 
similar susceptibilities to those we used, and b) the HRG 
approach - without any free parameter - is in very good 
agreement with the lattice calculation for T < 150 MeV. 
It is therefore unlikely that uncertainties in the suscepti¬ 
bilities will lead to dramatic changes in the sterile neu¬ 
trino PSDs. 


The validity of our neutrino opacities is much less clear 
- we have attempted to calculate them in as much detail 
as possible, but the hadronic parts still retain significant 
uncertainties due to the quark-hadron transition. We 
expect that opacities at high and low temperatures are 
well described by the rates of reactions involving free- 


quarks, and the lightest pseudoscalar and vector mesons, 
which are shown in Sec. IV For temperatures near Tqcd , 
we have considered two interpolation schemes (shown in 
Fig. 1^ that we expect might bracket the range of pos¬ 
sibilities. We have computed the sterile neutrino PSDs 
for both cases and shown that they are fairly robust to 
the choice used, especially for models with larger values 
of the average momentum (p/T). We leave the calcula¬ 
tion of self-consistent opacities through the transition to 
future work. Yet another approximation we have made 
is that of equilibrium distributions for all active species, 
which has been studied in a different context in Ref. |125) . 
We expect this to be valid at the temperatures relevant 
to the models we study. 

To compute the Hubble expansion rate and time- 
temperature relation, we have used the plasma’s equa¬ 


tions of state provided in Ref. |104j . which are obtained 
by matching to the lattice QGD results of |126j . As the 
former’s authors point out, this result is still uncertain 
at temperatures close to the quark-hadron transition. It 
would interesting to update their result with the latest 
lattice QCD computations, which suggest a lower tran¬ 
sition temperature m- We expect the uncertainties 
associated with the plasma’s equation of state to be at 
most similar in magnitude to those coming from the neu¬ 
trino opacity [20] . 


Another simplification we adopted is the semi-classical 
Boltzmann equation, which greatly facilitates our study 
of the oscillation-driven production. As mentioned in 
Sec. |HB[ the most general analysis considers the evolu¬ 
tion of a two-state density matrix, rather than phase- 
space densities. The validity of the semi-classical ap¬ 
proach rests on the assumption that collisions domi¬ 
nate the off-diagonal element of the Hamiltonian that 
is responsible for vacuum oscillations jSSHnZj. For 
typical momenta at the temperatures of interest, the 
ratio of these terms is A(T) sin^ 20/I?(T) ~ 0.6 x 
(T/lOO MeV)-®(ms/7keV)2(sin2 26 »/ 10 -“)i /2 ^he pro¬ 
duction of sterile neutrinos happens at temperatures 
above, but close to where these terms become compa¬ 
rable (note the ratio’s steep temperature dependence). 
Thus, we expect that the results in this paper are rela¬ 
tively unaffected by this approximation, but further work 
in this direction can settle this question. 


Finally, we examine the assumptions underlying the 
model itself, which were enumerated in Section |H A| If 
there is indeed an extra neutrino that is an electroweak 
singlet, it is not restricted to mix with only one fla¬ 
vor. However, the general case where the sterile neutrino 
mixes with all flavors introduces extra mixing angles, 
which cannot be constrained as easily from observations. 
The same can be said about the assumption of a lepton 
asymmetry in a single flavor. We briefly remark on the 
possibility of the sterile neutrinos mixing with electron 
or tau flavors instead. The redistribution of Sec. nni is 
almost identical for the cases with electronic and muonic 
lepton asymmetries, but is different in the tauonic case. 
This is due to the significantly larger mass of the corre- 
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spending charged lepton {rrir = 1.77 GeV m), which an¬ 
nihilates away at higher temperatures. Thus most of an 
input tau asymmetry ends up in the tau neutrino below 
T < 400 MeV, and the quark hadron transition does not 
impact the redistribution. The electron and tau neutrino 
opacities are different from the muonic case, and so is the 
balance between the thermal and asymmetry potentials, 
which affects the resonant momenta and ultimately the 
final dark matter PSDs - we leave for future work the 
possibility of sterile neutrinos mixing with those flavors. 


Also worth considering is active-active neutrino mix¬ 
ing, which does not conserve asymmetries in the individ¬ 
ual flavors. This was studied in Ref. [90], which showed 
that such asymmetries are frozen in at the temperatures 
of interest. An interesting possibility is to revisit this 


study and use the redistributed asymmetries of Sec. Ill 


to calculate the active neutrino self energies at this epoch. 


In conclusion, we find remarkable that sterile neutrino 
models that are in agreement with the X-ray excess have 
transfer function shapes that can signihcantly impact 
structure formation on subgalactic scales. Fixing the lep- 
tonic asymmetry to produce the right DM relic density, 
the resonantly-produced sterile neutrino transfer func¬ 
tion goes from ‘warm’ to ‘cold’ as the mixing angle is 
increased from small to large values. This indicates that 


upcoming X-ray observations |128l I129j and ongoing ef¬ 
forts to study small-scale structure can together cover all 
of the allowed mixing angle parameter space, and conse¬ 
quently confirm or disfavor the model. 
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Appendix A: Functions for high-temperature QCD pressure 


The functions defined in section III B I are |104j 


^2(J/,A) 

F4{y,fi) 


247r^ 


dx t 


x + y 


1 

pOO 

/ dx 

Jo 

X 

Stt^ 

_x + y 


oo ^ 

dx — 

X 

Jo 

X 

x + y\ 


1 


dxi / dx 


[fipiy/x + y - A) + npiy/x + y + A)] 

[fipW^ + y - A) + fipiVx^ + A)]) 

[hF(Vx+T/-A) + ^F(\AcTy + A)]: 

1 


(47r)4 ^xi + y^Jx 2 + y 

X I [fi^Wxi +y - {j)n-p{^/x 2 +y + A) + fi-p^xi + y + fl)nFiVx 2 + y - A)] 

X In 


W^i + yVx 2 + y + y- y/xiX 2 ' 

+ y + ^/x\X2_ 

[fiFiy/xi + y - A)«f(vAe2~+1/ - A) + fiFWxi -f y + jl)nFWx 2 + y + A)] 
'y/xi + yy/x 2 + y -y + VYTxi" 


X In 


'Jxi + y^x 2 + y -y- ^JxiX 2 


(Al) 

(A2) 

(A3) 


(A4) 


where nF(a;) = l/[expa; -|- 1] is the Eermi-Dirac distribution. 
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TABLE IV: Representative examples of two-particle to two-particle reactions that contribute to the neutrino opacity, along 
with their matrix elements written using a four-fermion vertex for the weak interaction (for reactions involving leptons and 
quarks) and three quark chiral perturbation theory (for reactions involving the pseudoscalar meson octet). For a reaction 
A ^ B -\- C , the momenta pi : i = 1, 2,3,4 within matrix elements are mapped to A, B and C respectively. 


Reaction 

SEI-Mp 

Remarks 

+ T~ ^ fl~ + Ot 

Reactions involving leptons 

128G| {pi ■P 2 ) {P3 - Pi) 


t'M + dw -f d 

Reactions involving quarks 

384 (pi • Pi) (p2 ■ p-i) 

T > Tqcd 

Ufj. + ^ TV'^ + 71 ° 

77^ -f TT’*' —>• -1- TV'^ 

Reactions involving the pseudoscalar meson octet 

8\Vud\^Gp [2p2 ■ {pi - P3)pi ■ {pi - ps) - (pi • P 2 ) {pi - Ps)^] 

4 (1 - 2sin^ G| [2p3 • (p4 -|-p2)pi • (p4 +P2) - (pi • Pa) (P4 +P2)^] 

T < Tqcd , s-channel ^ 
T < Tqcd, t-channel ^ 


Input neutrinos can produce quarks at low temperatures, T < Tqcd> for large CM energies in the s-channel. 

Input neutrinos can probe the quark content of mesons at low temperatures, T < Tqqjj, for large momentum transfers in the t-channel. 


Appendix B: Neutrino opacities: matrix elements and collision integrals 


In this appendix, we expand on the details underlying the neutrino opacities that were presented in Sec. EYl Under 
the set of assumptions presented therein, we add contributions from a large number of reaction rates, which we present 
in an organized manner in the rest of this section. 


1. Rates for neutrinos to go to two-particle final states 

We compute reaction rates for momenta and temperatures where we can integrate the weak gauge bosons out and 
approximate the weak interaction by a four-particle vertex. For tree level processes under this approximation, if one 
of the ingoing particles is a neutrino, one of the other particles is either a neutrino or a charged lepton belonging to 
the same generation. We classify reactions as leptonic or hadronic based on the nature of the remaining two particles. 

We now describe our calculation of these reactions’ matrix elements, and the associated contributions to the neutrino 
opacity. 


a. Matrix elements for two-particle to two-particle reactions 


It is a lengthy, but straightforward, task to enumerate all leptonic reactions that contribute to the neutrino opacity. 
Ref. mu lists a complete set of reactions at temperatures of a few MeV. Our calculations extend to higher tempera¬ 
tures, hence we also include reactions involving tau leptons. The reactions are enumerated in Table It is harder to 
study hadronic reactions in a consistent manner through the quark-hadron transition temperature, Tqcd- We adopt 
assumption ^ 3 of Sec. IV we neglect the strong coupling constant and its running at temperatures T > Tqcp, 


and hence calculate opacities with free quarks. We enumerate all reactions involving quarks in the same manner as 
the leptonic ones, but with standard model quark currents that couple to and W^. These reactions are listed in 
the upper half of Table |ll] Table |IV| shows squared and spin-summed matrix elements for a representative reaction 
involving only leptons, and for one involving quarks. 

The physical rates for hadronic reactions diverge from our calculated ones close to the transition, since the strong 
coupling constant is non-zero. Treating this self-consistently is beyond the scope of this paper. In the main body, we 
present results for a few unphysical interpolations through the transition. At even lower temperatures, T < Tqcd, 
we cannot use the free quark approximation. The most important hadronic degrees of freedom are the pseudoscalar 
meson octet, which are pseudo-Goldstone bosons associated with the spontaneous breaking of the axial part of an 
approximate SU (3) l x SU {2>)r flavor symmetry [130]. We use three quark chiral perturbation theory (3xPT) to write 
down the currents that couple to and W^, and through them evaluate the mesonic contribution to the neutrino 
opacity. 

Consider a 3 x 3 unitary matrix, U{x), which represents low-lying hadronic excitations at temperatures T < Tqcd- 
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We express U{x) in terms of the pion fields as follows: 

.7r“(a;)T“ 


U{x) = exp 
27r“(a;)T‘ 


2 i 


U 


1 

u 


(’ 


V 


u 

0 _|- 

V 2 K- 


a e [1,8], 

V2K+\ 


'/ 2 tt 


-TfO + ^r; 


VMo 


(Bl) 


(B2) 


71V 


where is an energy-scale associated with the breaking of the SU{3)a symmetry, and T“ are generators of SU{3). 
The most massive member of this octet, the rj meson, has a mass of to,; = 547.8 MeV [1]. We only use this prescription 
at T < 150 MeV, so these low lying excitations are sufficient to describe all relevant incoming hadronic degrees of 
freedom. 

In the framework of SyPT, the dynamics of the pion fields are described by an effective Lagrangian for U{x) coupled 
to matrix valued SU{3)l and SU{3)ii gauge fields and respectively. The first approximation to the Lagrangian 
is the lowest term in a derivative expansion: 


C = -^f^TT[D>^U^D^U], with D^U = df,U -il^U + iUr^. 
The gauge fields lfj_ and are Hermitian matrices, which we decompose as 

{l/r)^ = {l/rr.T‘^ + T A^)I, a e [1,8]. 


(B3) 


(B4) 


The fields Vfj, and Afj, are vector and axial-vector parts of Ifj, and r^. We identify the electroweak gauge bosons of 
the standard model, and [or equivalently, the underlying SU{2) x U{1) gauge fields and with 

elements of {A IVIA)^ by equating their action on the pion fields tt‘^{x) or the excitation U{x) via the right-hand 
side of Eq. (|B3|). The results of this procedure are 


92 AI = ll, 

92 Al = ll, 

92AI = ll + ^ll - + ^A^, 


1 


1 


exl„ — T r 4—^7 4— 

^ ^ ^ s/3 ^ y/3 ^ 


(B5a) 

(B5b) 

(B5c) 

(B5d) 


Here 52 is the SU{2) coupling constant, and e = 52 sin The overlap of with the anomalous A^ is an artifact 
of the omission of the charm quark from our analysis. If the charm quark were included (as it has to be to make the 
SM non-anomalous) then there is no coupling to an anomalous current. As long as there are no charm quarks (always 
true at the temperatures where 3xPT is applied) this current looks like that associated with the anomalous U{1)a- 
We read off the currents that couple to the gauge fields using the definition = dL j d{l°' jr'^ jV j A) and 

the Lagrangian of Eq. (B3). We transform the resultings to obtain the currents that couple to the SM electroweak 
gauge helds, and A^-. 
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(B6a) 

(B6b) 

(B6c) 

(B6d) 

(B6e) 


These currents agree with the leading order parts of the functionals computed in Ref. m- The lower half of Table 
[II| enumerates the allowed reactions involving pseudoscalar mesons. Table jlVj shows the squared and spin-summed 
matrix elements for two such reactions, computed using the currents listed in Eq. (B6) 
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A final complication is that S^PT, and the currents derived from it, are valid only when the momentum in the 
intermediate weak gauge bosons is low compared to the energy scale 47r/,r ~ 1 GeV [130) 1. The physical currents that 
couple to the SM electroweak gauge fields are continuous functions of this momentum; they approach the SM free 
quark currents for large momentum values. This manifests as the production of quarks in the large CM energy limit in 
s-channel reactions, and as ‘deep-inelastic scattering’ off the mesons’ quark content in the large momentum-transfer 
limit in t-channel reactions. These limits are important to consider at the higher energies for which we calculate 
neutrino opacities using SyPT (the total energy range is shown in Figures 8a and|^. 

We do not self-consistently compute these corrections to the currents, as it is beyond the scope of this paper. 
Instead, we modify the s-channel reaction rates in a phenomenological manner: we apply a cutoff in the CM energy 
at 1 GeV with a width of 50 MeV, below which we use the 3xPT currents, and above which we use the SM free quark 
currents. We do not incorporate any corrections to t-channel reactions; this would involve some knowledge of the 
parton distribution functions for the mesons involved. 

We observe that the squared and spin-summed matrix elements for tree level processes, be it for processes involving 
leptons and free quarks (computed using the SM currents), or for those involving pseudoscalar mesons (computed 
using SyPT), are at-most quadratic functions of the Mandelstam variables. This greatly facilitates a semi-automated 
computation of the two-particle to two-particle reactions’ contribution to the neutrino opacity, which we very briefly 
describe next. 


b. Rates for two-particle to two-particle reactions 


Consider a general two-particle to two-particle reaction, j / q , -|-A — >■ B -\-C, that consumes a massless input neutrino, 
Ua- The particles A^B and C can all be fermions (leptons or quarks), or contain a pair of bosons (pseudoscalar 
mesons). We expand Eq. 0 to write down the following expression for the scattering rate as a collision integral: 

^ ydWpBrf3pc(27r)4<5(p,^ + PA- PB - Pc)SY,\M\‘^ fA{EA)ilT fB{EB)){lT fc{Ec)), (B7) 

where the symbol d^p is shorthand for the Lorentz invariant phase space volume element d^p/[(27r)^2E(p)], the 
symbol ^ |A4p is the absolute value of the matrix element squared and summed over all spin states, S' is a symmetry 
factor for identical particles in the initial and/or final states, and the f{E)s are appropriate Bose-Einstein/Fermi- 
Dirac phase space distributions depending on the statistics of the particles, with plus and minus signs for bosons and 
fermions respectively. 


We follow the treatment in Ref. [Ill] to reduce the nine-dimensional phase space integral of Eq. (B7) to a numerically 
manageable three-dimensional integral over the variables |pa|, |pb| and pb =Pb ’Pua- This procedure involves using 
the delta function to perform the integral over pc, and using the form of the matrix elements for tree level processes 
to analytically perform the integral over ra = Pa ■ Pi^ ■ We refer the reader to Ref. m for more details. The form 
of the matrix elements also lends itself to easy parameterization in terms of a small number of classes; along with 
the procedure described above, this enables a simple numerical implementation of the calculation of these reactions’ 
contributions to the neutrino scattering rate. 


2. Rates for neutrinos to go to one-particle final states 

We must also consider the contribution to the neutrino interaction rate, from interactions with two-particle 

final states (“fusion” or inverse decay). A four-fermion interaction such as the weak interaction (at E <C mwTBiz) 
can produce such a final state in two ways. One, applicable at T < Tqcu, is two-body fusion to produce a meson, e.g. 

-I- —)■ 7r+. The other is the ‘three-body fusion’, e.g. -I- T-e -I- e“ —tp~. By construction, these fusion processes 

are the inverse of a decay process. We describe our treatment of these processes in the rest of this section. 


a. Kinematics of two-body fusion 

A two-body fusion process must involve a meson in either the initial or the final state, and - if it is to absorb a 
neutrino - must then be semi-leptonic. The neutral current processes of this form (e.g. a neutral meson is created by 
the fusion of BaPa t^'^) are helicity-forbidden and have zero amplitude at tree level. The charged current processes 
can have either the meson in the initial state and the charged lepton in the final state (e.g. K~Vt -a t~) ov the meson 
in the final state (e.g. -A 7 r+). The “charged lepton in the final state” case is possible only if the charged lepton 
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is more massive than the meson, i.e. if that lepton is a t; at T < Tqcd this not energetically feasible for typical 
values of the incoming particles’ momenta, since tUt Tqcb- Therefore, for the rest of this section, we focus on the 
problem of a charged meson in the final state. The reaction is 

i/a + a+ A+ (B8) 


where a = e oi fj, and A = tt or K. We are interested in the thermal absorption rate Tfugion for the neutrinos as a 
function of temperature T and neutrino energy Ei,. 

The simplest solution to this problem is to calculate the rate of the inverse reaction of Eq. (B8) and use detailed 
balance. In thermal equilibrium, there is a rate of decays given by 


dN 
dV dt 




47rp^ dp A 
(27r)3 


dp! niA 


fA(EA)[l-UjE,J][l-UE^)], 


(B9) 


where p' is the cosine of the angle of emission of the neutrino in the rest frame of the A+, and the factor of vriAlE a 
is the inverse-Lorentz factor that accounts for the longer lab-frame lifetime of A'^ at high energies. The degeneracy 
factor is 5 a = 1 for pions and kaons, but we include it for later use with heavy mesons. 

With the help of relativistic kinematics, we see that the lab-frame neutrino energy is 


= ^ ( 1 - 


TOt 


{Ea +PA^J-'), 


(BIO) 


so that Eq. (B91 can be re-written in terms of a rate of decays per unit volume per unit neutrino energy: 

47rp^ dp A 1 


dN 

dV dt dE^ 




ruA 


(27r)3 (1 - m1lm\)pA Ea 


The kinematically allowed range of A~^ momenta is given by 


or - using pA = E\ — m\ and with some algebraic manipulation - 

(1 - rn^/m\)‘^m\ + AE^^ 


/a(Ea)[1 - UAEuJ][l - /a(E„)]. (Bll) 
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Turning the integral into one over the energy Ea of A'^ gives 
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Now this should equal the fusion rate of neutrinos, which is 

dN AnEf 


dV dt dE^ 


(27r)3 
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We therefore conclude that 


rfusion ) 


(1 - mllm\)El^ 


r fA{EA)^-^^f^^[l-UEa.)]dEA. 
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(B16) 


Next we substitute in the Bose-Einstein or Fermi-Dirac phase space distributions, and note that Ea = Ea — E ^^, 
yielding 


Tfusion {E.^) 
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Here we have defined the dimensionless integral 
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Recall that i?_A,min is a function of and is given by Eq. (B13). We then achieve the final simplification: 

J- fusion j - /- I-. /rr.v in- 


;(1 + e -®- c / r )£;2 


1 _ e-(’^^™i+4-E2^)/(4n_E,„T) 


(B18) 


(B19) 


where r; = 1 — Note that the numerical calculation of the logarithm must be treated carefully since for 

Eyi^min ~ Ei,^ ^ T we are taking the logarithm of a number that is very close to 1. For calculational purposes, we 
replace the logarithm in Eq. (B19) by a truncation of its Taylor expansion at the fifth order wherever the argument 
deviates from unity by less than e = 10“^. 


b. Rates for two-body fusion processes 


The rate parameters for the key two-body fusion reactions are shown in Table 

Some parameters were not available in the Review of Particle Properties. Key among these are the decay parameters 
for the weak decays of the p, a;(782), and K*{892) vector mesons. These mesons are actually extremely broad 
resonances, and their principal decay mode is into lighter mesons. Electromagnetic and especially weak decay modes 
are less common. Note that the helicity suppression arguments that forbid e.g. tt^ —>■ i^eZ/e do not apply to the vector 
mesons. 

The relevant decays of the charged vector mesons p~^ and Rr*(892)+ can be obtained by noting that the virtual-IT 
diagram results in an effective vertex 


/lefF 3 ^ 
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2-\/2 sin 0VV 




1 e 

‘2\p2 sin0u^ 


yud‘rnpfp\e{p )]^-|-h.c.. 


(B20) 


where fp is the p decay constant, and \e{p~y\p is the polarization of the “on-shell” p meson. The advantage of this 
Lagrangian is that the well-measured decay —)■ p'^Dr is related to the decay of to a charged lepton and a 

neutrino (predicted branching fraction ~ 2 x 10“^^). From a tree-level calculation with this Lagrangian, we infer a 
ratio 


r(p+ ^ p+iyp) _ 2 mUml - w^)^( 2 mg -b w^) 
r(r+—)• p+Pt) — mp)2(2TO2-I-m^) 


(B21) 


These rates are included in Table IVl 
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TABLE V: The parameters for reactions that go into Eq. ( |B19 1. Reactions relevant for the neutrino opacity are shown; 
antineutrinos are similar. Particle masses are obtained from the Particle Data Group. Decay partial widths are obtained from 
the sources indicated. All reactions in which a neutrino can produce a hadronic resonance below 1 GeV are included. 


Reaction 

ruA 

QA 

V r™'" 

reverse 

Rate method 


MeV 


MeV 





Reactions involving the pseudoscalar meson octet 

Vs + ^ 

139.57 

1 

0.999987 3.110 x 10"^® 

PDG 

Vp+ p+ ^ 7r+ 

139.57 

1 

0.4269 2.528 x 10"^'‘ 

PDG 

Vs + e+ ^ K+ 

493.68 

1 

0.9999989 8.41 x 10"^® 

PDG 

f^fJ. + 

493.68 

1 

0.95419 3.38 x lO"^'^ 

PDG 




Reactions involving vector mesons with nonzero isospin 

Vx + Vx ^ p° 

775.26 

3 

00 

X 

0 

1 

to 

Average of r decay and e+e“ — >■ p°; assumed isospin 5'17(2|j 

Vs-\- e'^ ^ 

775.26 

3 

0.9999996 7.00 x 10"“ 

r decay 

Vp+ p'^ ^ 

775.26 

3 

0.98143 6.80 x 10"“ 

T decay 

Vs + e+ ^ if* (892)+ 

891.66 

3 

0.9999997 5.45 x 10"^^ 

r decay 

vp + p+ ^ K*{8Q2)+ 

891.66 

3 

0.98596 5.33 x 10"^^ 

T decay 




Reactions involving vector 

mesons with zero isospin 

vx -\- vx ^ u>{782) 

782.65 

3 

1 7 x 10"^® 

e+e“ -A- lo{7S2)-, assumed quark content {uu + dd)/y/2 


^ The r decay gives 1.01 X 10 and the e+e —> computation gives 9.5 X 10 


The rate for > vx^x can be obtained by replacing the terms in Eq. (B20) with the Z couplings and propagator: 




4 sin 9w cos 9w 


- l^)i’x 
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sin 9\^ cos 9w 
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(B22) 


Here the coupling of the vector Z current to the meson was related to the coupling of the vector W current to the 
p'^ using isospin SU{2) symmetry. The conclusion is that 


8m^ 


m„ 


r(p° -» VeVe) _ _ 

T{t+ ^ p+Dr) 3(to 2 — TOp)2(2m2 + m^) 


(1/2 — sin^ 9wY 


IK 
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(B23) 


and similarly for the other two flavors. In this relation, we have used that mw = mz cos 9\y- 
The same decay rate can be obtained by taking the ratio 


r(p0 




TipO 


+ #3—1 


e'^e 
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2 sin 6\y cos 9^ 



1 {0\{9u,vuj'-u + gd,vd'y’-d + gs,vsj''s)\a° 

sin0vvcos0w (OKluy®!!— hd^d— is7®s)|/9°) 


Note that the reaction 


with guy = 4 ~ 5 sin 9w 


(B24) 

■ e~ is mediated principally by the photon instead of the Z. The ^ is due to the fact 
that the photon couples to right-handed as well as left-handed electrons. The factor of rria/mp is a combination of 
kinematic factors appropriate if the decay constant is the same for all octet members. The factor of 1/(2 sin 9w cos 9w) 
is the ratio of the Z coupling to to the photon coupling to ez, (or e^). The factor of m}plm\ is the ratio of Z to 
photon propagators. The last term is the ratio of Z coupling to to 7 coupling to 
and gay = gsy = + 5 sin^ 9\y. For the the last term can be computed using isospin symmetry. 

The agreement between r decay and e“'"e“ —>■ is good: the former predicts a partial width for p^ —> vZVf. of 
1.01 X 10“^^ MeV, and the latter predicts 9.5 x 10“^^ MeV, a difference of only 6%. The average is shown in the 
table. 

No similar decay is allowed (i.e. it is not possible with a single intermediate vector boson propagator) for the 
AT* (892)° or .^*(892)° mesons because the current that couples to the Z cannot change strangeness. 

It is less clear how this procedure should be applied to the uj{782) meson, which has no isospin. One might 
approximate it as a pure {uu + dd)/\/2 state (i.e. with no strange quark), and repeat the argument used for the 


e^e 


calculation with Eq. (B24|). This result is shown in the table; it is much more uncertain than the calculation 


P _ _ 

for the p° since mixing with ss is allowed. Nevertheless, the small rate for a;(782) production (as compared with p°) 
suggests that it leads to an overall small correction to neutrino opacities. 
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c. Three-body fusion processes 

The final set of reactions that contribute to the neutrino opacity are three-body fusions. As earlier, these reactions 
can be either leptonic or hadronic in nature. We adopt the prescription outlined in Sec. B 1 a for the hadronic 
reactions. Given the hadronic and leptonic currents coupling to the SM electroweak gauge bosons, we can enumerate 
all three-body reactions that contribute to the neutrino opacity in same manner as earlier. We must keep in mind 
the kinematic constraint that the rest mass of the product must be greater than that of the reactants. 

The matrix element for any three-body fusion reaction is related to one for a two-particle to two-particle reaction 
by crossing symmetry. Thus, we do not need to compute any new matrix elements for this section. However, we need 
to modify the treatment of the kinematics from the previous case. Consider a general three-body fusion reaction, 
Va A B ^ C. The scattering rate for an input neutrino energy is given by the collision integral: 

J d^pA(fpBd^Pci2Tr)'^S{p^^ +PA+PB -pc)S ^ \M\^ fA{EA)fB{EB){^ T fc{Ec))- (B25) 


All the symbols are defined identically to Eq. (B7). The procedure to reduce the dimensionality of this integral is 
exactly analogous to that in Sec. Bib and Ref. Uni, with one important difference. The variables final y left to 
numerically integrate over are, as earlier, |pa|, |pb| and pB = Pj/„ ■ Ps- If we consider the integration domain for the 
two-particle to two-particle case, for a given value of |pyi|, energy constraints allow a maximum value of |p_b|. For a 
three-body fusion, |pb| has no upper bound, which greatly expands the allowed phase-space. With this caveat, the 
rest of the procedure proceeds as it did for the other case. 
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